Turbulent sensible and latent heat exchanges play an important role in melting snow covers, contributing 30–40% of overall melt energy with daily values reaching over 50% on warm, cloudy days (Morris, 1989). The spatial variability of these turbulent fluxes across a basin and the relative importance of the differences is not well known. This paper specifically addresses small-scale variabilities in sensible and latent energy fluxes related to topographically induced wind speed variations. A simple wind model was used to simulate topographic effects on the surface wind field. Hourly wind observations were areally distributed by the model and used to calculate spatially variable sensible and latent turbulent heat fluxes for a small (63 km2) research catchment dominated by open tundra vegetation. Simulations showed that, even though the study area is characterized by relatively low relief (average slope 3°), the small-scale sensible and latent heat fluxes varied considerably throughout the basin. The resulting variations in snowmelt rates play an important role in the development of a patchy snow cover. Overall, turbulent fluxes within the research area varied by as much as 20% from the mean, leading to differences in potential snowmelt of up to 70 mm snow water equivalent over the entire melt period.
Introduction
Arctic snowpacks typically become patchy early in the melt season, with snow covered areas and bare ground coexisting throughout the spring. This phenomenon can be attributed to complex interactions between a highly variable end-of-winter snow cover and spatially variable snowmelt energy fluxes. The patchiness of the decaying snowpack further enhances the small-scale melt variability. Determining the extent and location of snow covered and snow free areas is crucial for predicting the quantity and timing of meltwater release from different parts of a basin and the energy exchange between the land surface and the atmosphere.
While arctic end-of-winter snow covers are generally continuous, they are characterized by high snow water equivalent (SWE) variability due to redistribution of snow by blowing snow events. Source and sink areas of blowing snow are mainly determined by topography, vegetation, and the wind direction of each individual blowing snow event. Blowing snow events also cause a significant sublimation loss of snow. Liston and Sturm (2002) and Pomeroy et al. (1997) reported average sublimation losses over the winter period of 9–22% of annual snowfall during blowing snow events in arctic catchments, with much higher values being simulated on windward slopes.
Spatially variable melt energy fluxes are superimposed on this variable snow cover leading to the quick development of a patchy snow cover. Incoming solar radiation is spatially distributed due to topographic effects. Studies have shown the importance of these small-scale variabilities for snowmelt in mountainous (Marks et al., 1999) and low relief areas (Hinzman et al., 1992; Pohl and Marsh, 2006). Pohl et al. (2006) found that small-scale solar radiation variabilities may lead to differences in melt of up to 50 mm SWE within a small arctic basin.
Another major component of the snowmelt energy balance are turbulent fluxes of sensible and latent heat. Turbulent fluxes are especially important in open, wind-swept environments (Male and Granger, 1981). As the air temperature increases above 0°C, the temperature difference between air and snow, which has a constrained surface temperature, becomes positive and sensible heat fluxes are directed toward the snow cover, contributing to melt (Marks and Dozier, 1992). Hinzman et al. (1992) noted that while condensation is quite common over melting arctic snow covers, latent heat is mostly directed away from the snow surface. This evaporation leads to a significant loss of snow of up to 33% that can, therefore, not contribute to runoff (Eaton and Wendler, 1982). During spring melt periods, latent heat often mirrors sensible heat, reducing the amount of melt energy available from turbulent fluxes.
Turbulent fluxes further contribute to snowmelt by counteracting the nocturnal radiative energy loss of the snowpack, which is especially high in clear nights (Martin and Lejeune, 1998). The turbulent fluxes to or from the surface depend strongly on atmospheric conditions, and especially on the temperature profile. An atmosphere is considered neutral when turbulence is not affected by changes in air density caused by temperature gradients. Under low wind conditions, however, temperature gradients can develop and result in an unstable or stable atmosphere. Unstable atmospheric conditions with warm air underlying cold air enhance turbulence and therefore increase the turbulent energy exchange, while stable conditions with cold air underlying warm air suppress turbulence. Studies have shown that temperature inversions (stable atmospheric conditions) often develop during nighttime over melting snow covers, resulting in reduced turbulent fluxes especially early in the day (Martin and Lejeune, 1998). Overall, sensible heat input to the snowmelt energy balance can be up to twice the net radiation during cloudy, warm, and windy periods (Pohl and Marsh, 2006) but is very often around 40% of net radiation, while contributions of latent heat can be of the same order as net radiation but are typically around 10% (Morris, 1989).
Turbulent fluxes can be measured with eddy covariance techniques. However, due to the complexities of these techniques (Blanford and Gay, 1992), turbulent fluxes are commonly calculated using bulk aerodynamic formulae. This approach requires surface temperature information and observations of wind speed, air temperature, and relative humidity at at least one atmospheric level. Bulk transfer coefficients for momentum, heat, and water vapor are calculated, initially for a neutral atmosphere, and subsequently modified to account for stable or unstable atmospheric conditions (Moore, 1983).
The topographic influence on airflow near the surface has been widely recognized in the literature. However, there have been few studies addressing how these flow variations affect turbulent energy fluxes over melting snow. In fact, many spatially distributed snowmelt studies assume constant wind speeds over the entire model domain (Bloeschl et al., 1991; Hinzman et al., 1992) or distribute wind speed using linear elevation gradients (Susong et al., 1999). Topographic influences on the surface wind field can be simulated with a wide variety of models. The most complex simulations apply a physically based, fully atmospheric model or a boundary layer circulation model (e.g., Liston et al., 1993). Another approach is the use of an atmospheric model in which only mass continuity is satisfied (e.g., Ross et al., 1988). Computationally simpler methods use observed data to define empirical wind-topography relationships that are then employed to develop a spatially variable wind field from initial point measurements (Marsh, 1999; Liston and Sturm, 1998; Purves et al., 1998). Other studies rely on utilizing meteorological stations located in representative topographical locations, and relate the observations to similar locations within a study area (Woo and Young, 2004).
This paper presents a study of the small-scale variability in turbulent fluxes introduced by wind-topography interactions and their effect on snowmelt. It is the objective of this study to show the amount of variability present in the sensible and latent heat fluxes, even in the relatively gentle terrain of the study area, as a result of distributed surface winds. Furthermore, interactions between wind dependent spatial distribution of the end-of-winter snow cover and spatially variable turbulence fluxes are shown and their role in the formation and disappearance of a patchy snow cover is discussed.
Study Area
The study was conducted in Trail Valley Creek (TVC), a National Water Research Institute research basin in the Northwest Territories of Canada, located at 68°45′N, 133°30′W. The area is characterized by gently rolling hills with a small percentage of deeply incised river valleys. Elevation ranges from 40 to 187 m a.s.l. (Fig. 1) with an average elevation of 99 m a.s.l., a mean slope of 3°, and a maximum gradient reaching 33°. The region lies at the northern edge of the forest-tundra transition zone, with tundra vegetation dominating much of the upland areas, and some shrub tundra and sparse black spruce (Picea mariana) forest on hillslopes and in the valley bottoms (Neumann and Marsh, 1998).
The prevailing winter and spring winds in this area are from the northwest to northeast. The lack of winter melt/freeze cycles, in conjunction with frequent high wind speed events, leads to extensive snow redistribution due to blowing snow. Detailed studies (Essery et al., 1999; Pomeroy et al., 1997) have shown that snow is removed mainly from open tundra areas, and is deposited in regions with taller shrub and forest vegetation. Essery et al. (1999) showed that snow erosion occurs preferentially on windward slopes, while accumulation in snow drifts is greatest on leeward slopes relative to the dominant winter wind directions. Overall, drift areas cover about 8% of TVC (Marsh and Pomeroy, 1996). Snow redistribution leads to a highly variable end-of-winter snowpack in TVC. For example, in 1993 SWE expressed as a percentage of annual measured snowfall varied from 54% for tundra areas to 180% for areas of taller shrub and forest vegetation, to 419% for drifts (Pomeroy et al., 1997).
Snowmelt in the study region typically occurs in May and June. The present study was conducted from early March to early June 1998 and 1999 as part of the Mackenzie GEWEX (Global Energy and Water Experiment) Study (MAGS)/Canadian GEWEX Enhanced Study (CAGES) program. Two permanent meteorological stations were situated in the basin (Trail Valley Main Met [TMM], and Trail Valley Upper Plateau [TUP]). Wind speed was also observed at two additional temporary met stations within the basin (Trail Valley Deep Snow [TDS], and Trail Valley Bare Ground [TBG]) (Fig. 1). A digital elevation model (DEM) of a 14 × 12 km area that includes TVC was used to model the topographic effects on surface winds. The DEM, consisting of 105,000 grid cells with a grid increment of 40 × 40 m, was obtained by digitizing 1:50000 National Topographic Survey maps.
Methods
MODEL DESCRIPTION
Two models were considered in the present study to simulate topographic influences on the surface wind field. The MS3DJH/3R model developed by Walmsley and Taylor (Walmsley et al., 1982, 1986; Taylor et al., 1983) represents a moderately complex windflow model. It is a linear model developed from theoretical analysis of neutrally stratified boundary-layer flow over low hills. This model was used by Essery et al. (1999) to drive a distributed blowing snow model for TVC. Essery (2001) noted that despite its linearity, the model is still computationally demanding, especially when applied to large areas at high resolutions. A simpler approach was introduced by Liston and Sturm (1998). This “Liston model” uses point wind speed and wind direction observations in conjunction with empirical wind-topography relationships. Observed wind speeds and directions in the model domain are interpolated to the model grid and subsequently multiplied with empirically based wind weighting factors (WWF) to account for topographic variations. The model also includes an algorithm (Ryan, 1977) to modify the initial wind direction. The algorithm uses the initial wind direction, terrain slope, and azimuth to calculate a diverting factor for each grid cell of the DEM. This diverting factor is then added to the initial wind direction to obtain a terrain-modified wind direction (Liston and Elder, 2005). However, due to the relatively gentle nature of the study area, this yielded a maximum deviation of less than 6°. Thus, we assume wind direction to be unaffected by topography for the remainder of this paper.
In the present study, we use results from Essery et al. (1999), based on the Walmsley-Taylor model, to set the empirical parameters needed for the Liston model (Liston and Sturm, 1998), which is subsequently used to derive the WWF needed to drive the energy balance calculations. The Liston model is used to optimize computational efficiency over our high resolution domain, while retaining as much information on spatial variability as possible. In addition, the simplicity of the Liston model in representing a spatially distributed surface wind field makes it easily applicable to other sites for use in a variety of modeling studies. The model uses a DEM to calculate topographic slope and curvature factors for each grid cell with respect to an observed wind direction. The model reduces wind speeds on leeward slopes, while accelerating winds on windward slopes. It further assumes that convex slopes experience higher wind speeds, while wind speeds are decreased over concave slopes. Accordingly, the lowest wind speeds are found on leeward concave slopes, while the highest wind speeds are simulated on steep convex windward slopes. WWFs are determined for each individual grid cell, and are dependent on wind direction from the computed slope and curvature factors (Ωs, Ωc) and additional empirical weighting factors (γs, γc) from the equation (Liston and Sturm, 1998)
The empirical factors (γs = 1 and γc = 3) were chosen for the present study such that (when multiplied with Ωs and Ωc) curvature and slope exert approximately the same influence on WWFs and that WWFs close to ones simulated for the area by the MS3DJH/3R model were computed. Figure 2 shows the WWFs for TVC for a wind direction of 45°. A distributed wind field is subsequently determined by multiplying an initially uniform, observed point wind speed from a met station within the model domain (at a location with WWF of 1.0) by the computed WWF for the respective observed wind direction.
TURBULENT FLUX CALCULATIONS
The turbulent flux simulation was run for 33 days during the spring of 1998 and 1999. A bulk transfer method (Dunne et al., 1976; Moore, 1983) was used to estimate sensible (Qh) and latent heat (Qe) fluxes over snow covered areas from the following equations:
where ρ is air density in kg m−3, cp is the heat capacity of air in J kg−1 °C−1, U is wind speed in m s−1, Ta is air temperature in °C, Ts is snow surface temperature in °C, Lv is latent heat of vaporization in J kg−1, Mv Ma−1 is the ratio between water-vapor and dry-air molecular weights, P is atmospheric pressure in Pa, and ea and es are the water-vapor pressures of the air and the surface in Pa, respectively (Martin and Lejeune, 1998). Dh and De are the dimensionless bulk transfer coefficients which depend on instrument height and surface roughness and are calculated after where k is von Karman's constant, z is instrument height in m, and z0 is surface roughness height in m.
The same method is used in many land surface schemes of atmospheric models (e.g., CLASS) and snowmelt energy models (e.g., SNTHERM) (Jordan, 1991; Verseghy et al., 1993). Air temperature and humidity were measured at TMM and TUP at 2.5 m above the ground surface (1.1 m in 1998). Instrument heights were adjusted daily as the snow cover underneath the instruments melted. A surface roughness of 0.002 m for TVC was estimated from wind profiles (Marsh and Pomeroy, 1996). This value is in good agreement with others reported for melting snow covers in flat terrain (Moore, 1983; Morris, 1989). Transfer coefficients were corrected for stable and unstable atmospheric conditions after
where Qhs is the stable bulk transfer coefficient, Qhu is the unstable bulk transfer coefficient, and Ri is the Richardson number used as indicator of atmospheric stability. Ri is computed after where g is the gravitational constant (9.81 m s−2), and z0 is the surface roughness in m. The same correction procedure was applied to De.
Equation 2 is often used for melt conditions with Ts set to 0°C. However, several studies (Marsh and Pomeroy, 1996; Martin and Lejeune, 1998; Luce et al., 1998) have demonstrated that the use of a variable snow surface temperature greatly improves the accuracy of the computed turbulent fluxes as it properly accounts for periods when surface temperature falls below 0°C. Hourly snow surface temperatures were estimated by assuming that snow surface temperature lagged air temperature and including a calculation of radiative cooling (Jordan, 1991; Marsh and Pomeroy, 1996). This procedure was shown to work well when compared to measured surface temperature values (infrared probe) at the TDS location within TVC (Fig 3). Wind speed was measured at the same height as air temperature and humidity, and then distributed over the entire model area using the windflow model. The resulting spatially variable wind field was used to calculate hourly sensible and latent heat fluxes for every grid cell of the DEM. Positive fluxes indicate an energy gain by the snow surface.
Wind speed was the only variable used for the calculation of the turbulent fluxes that was spatially distributed for this study. There are three other climatological variables that affect turbulent fluxes: barometric air pressure, humidity, and air temperature. Over an area as small as our study basin (63 km2), there is no reason to believe that barometric air pressure varies significantly. The same is true for humidity. Relative humidity (RH) is measured within the basin at two stations TMM and TUP located in the northeastern and northwestern corners of the basin (see Fig. 1). Over the modeling period, the difference in averaged measured RH at the two stations is 2%, which is less than the measurement errors of the instruments. The error in our turbulent flux calculations was further minimized by averaging the values of the two stations. Air temperature over a small area generally varies mainly with elevation. Only very little variation is expected for a basin with a low elevation range such as TVC (40–187 m). In the spring of 1999, air temperature was measured in the basin at two permanent locations (TMM and TUP) and two temporary stations (TBG and TDS). Again, an average of the two permanent stations was used in the turbulent flux calculations. These two stations cover a large portion of the overall elevation range of TVC being located at 70 m (TMM) and 170 m (TUP). Over the model period, the temperatures averaged 0.20°C for TMM and −0.33°C for TUP, with −0.07°C being the average of the air temperatures used in the model calculations of the study. This means that the averaged air temperature used for the calculations differed by only 0.27°C from the two meteorological stations, which is lower than the typical instrument inaccuracies of ±0.2°C. The comparison of all four stations shows that average air temperature differences in TVC are less than 0.6°C. An analysis of the variability resulting from these temperature changes shows that they are much smaller than the ones caused by the distributed wind fields outlined in the following section.
Results
FIELD OBSERVATIONS
To test the wind models, a field experiment was conducted in TVC to study the airflow around a low hill during the spring and summer of 2002. Six cup-anemometers with wind vanes were placed on a hill in the northwestern part of the basin (Fig. 1). Figure 4 shows the anemometer locations and the WWFs for the hillslope study site for a wind direction of 45°. Table 1 gives additional information on the placement of the anemometers. The data, averaged over 15-min intervals, was divided into a mostly snow covered (1 May to 1 June) and a mostly snow free period (2 June to 7 Sept.), in order to separate the data by surface roughness (i.e., snow/no snow). The snow cover during the May period leads to a relatively uniform surface roughness at all locations. The summer ground surface at the study site is characterized by low tundra vegetation except around station 5 which was located near some low shrubs. The locations were chosen to represent a wide range of topographic settings including: station 1 on the hilltop; station 2 on the upper reaches of the north-facing hillslope; stations 3, 4, and 6 at the bottom of the hillslopes facing north, east, and west, respectively, with 3 in an east-west–running valley, while 4 and 6 had open flat fetches to the east and west; and station 5 on the middle part of the south-facing slope (Fig. 4).
The observed wind data were grouped into the eight principal wind directions (N, NE, E, etc.). To study and compare wind speeds at each location, wind speed ratios (WSR = wind speed at station X/wind speed at station 1) were calculated for each site and each wind direction. The WSR for all anemometers and each wind direction for spring and summer periods are shown in Table 2.
On average, the highest wind speeds were observed on the hilltop for all wind directions. The results clearly show a pronounced windward-leeward effect around the hill. Wind speeds at station 2 were almost equal to the hillcrest location for northerly and northeasterly winds while being 10% lower for southerly and southwesterly winds. This shows that even locations in close proximity to the hilltop and with a relatively small elevation drop exhibit clear windward-leeward differences. This effect is even more pronounced for station 5 located in the middle of the south-facing hillslope. Here WSR for southerly to southeasterly winds are close to 1.0, while northerly to northwesterly winds had WSR of approximately 0.8. Wind speeds around the bottom of the hill are generally much lower than at the top while still showing strong windward-leeward effects (see east winds vs. west winds at station 6). The results for station 3 demonstrate the funneling of the wind in valleys, with highest WWF being observed for easterly and westerly winds parallel to the valley orientation.
Overall, the experiment shows that a relatively small hill (maximum elevation difference 54 m) modifies the surface wind fields considerably, with maximum wind speed differences between 24% and 39%, depending on wind direction.
MODEL VALIDATION
The Liston model was run on an hourly basis for the period from 5 May to 6 June, 1998 and 1999 (JD 125–JD 157), with wind speed and wind direction data from the TMM met station. TMM is located in an extensive flat area and therefore has a WWF very close to 1.0. WWFs were determined for each grid cell of the DEM for the eight principal wind directions at 45° increments from 360° (north) to 315° and observed wind directions (from TMM) were rounded to the nearest of these directions. The model results were then compared to observed wind speeds for the grid cells in which the three additional wind measurement sites (TUP, TDS, and TBG) were located (Fig. 1). Figure 5a shows the comparison for the TUP site. Overall, modeled wind speeds at all three sites compared fairly well against observed values, with R2 values ranging from 0.71 to 0.76. Unfortunately, this comparison did not provide a good test of the model since all four stations were situated on relatively flat terrain with a narrow range of WWF, from 1.02 to 0.94. The analysis, however, proved that flat areas throughout the basin have similar wind speeds regardless of elevation differences.
To better validate the Liston model, model predictions were compared to observational data from the 2002 field experiment. These field data covered a range of WWF from 1.10 to 0.83. Since wind observations at TMM were done at different heights above the ground surface and with different anemometers than those used in the field experiment, it was decided to compare the observations of the field stations against each other rather than against TMM data. Wind speeds and directions at station 1 (hilltop) were used as measured values (i.e., WWF = 1.0), and modeled WWF (see Table 2) were used to simulate expected wind speeds for the locations of stations 2 through 6. Figure 5b shows the scatter plot of observed versus simulated wind speeds for the spring period for station 5. Again, modeled values for all stations agreed well with measured data, with R2 values ranging from 0.82 to 0.97.
The experiment also illustrates areas of the model that could be improved. While areas of lowest wind speeds were correctly predicted to be on the bottom (concave) parts of the leeward hillslopes, the model did not capture the areas of highest wind speeds accurately. It assumes them to be on the steepest convex parts of the windward hillslopes, while the experiment suggests that the actual hilltop is the windiest location. The model does predict higher than average wind speeds at the crest of hills and ridges due to the convex curvature of these features. The curvature routine, however, does not accelerate the winds at those locations enough. A simple solution would be to use a higher empirical curvature weighting factor γc. This, however, would lead to unrealistically low wind speeds in the concave narrow valleys of TVC. To fix the problem, it might be necessary to compute curvature at a larger scale or to calculate the average curvature over a larger area around the point of interest. Another possibility is the use of separate γc for convex and concave areas.
TURBULENT FLUXES
As mentioned, the model was run on an hourly basis for the time period between 5 May and 6 June, 1998 and 1999. Figures 6a, b, c show the meteorological conditions for the modeling periods. The spring of 1998 was dominated by warm air masses leading to a large turbulent heat flux contribution to snowmelt, with the sensible heat flux accounting for 148 MJ m−2, while condensation added another 31 MJ m−2 (Fig. 7a). In 1999, daytime temperatures rose above the freezing mark after May 5, leading to positive sensible heat fluxes. A cold spell from May 15 to May 21 reduced the amount of sensible heat being transferred to the snowpack before temperatures and sensible heat fluxes rebounded. The latent energy exchange was dominated by sublimation during the early part of the melt. Significant contributions to snowmelt energy by condensation were only evident after May 29 (Fig. 7b). Over the whole period, sensible heat flux delivered 72.3 MJ m−2 to the snowpack, while the pack lost energy (4.2 MJ m−2) due to latent heat flux. Average wind directions (Fig. 8) and wind speeds (within 7%) for the two years were nearly identical, with northerly to northeasterly winds dominating both modeling periods.
The calculation of latent heat fluxes is fairly sensitive to measured values of relative humidity. A 5% change in relative humidity would result in an uncertainty in computed latent heat fluxes of ±8 MJ m−2 (or an hourly average of ±2.8 W m−2) over the study period. These uncertainties, however, are small compared to sensible heat flux values. Furthermore, the comparison of observed humidity values at TMM to sling psychrometer measurements showed a fairly good agreement (R = 0.91).
SPATIAL VARIABILITY IN THE TURBULENT FLUXES
Figure 9 shows the spatially distributed turbulent heat flux for 28 May (JD 148) at 19:00 h. Air temperature was 9°C, RH was 67%, and wind speed (for WWF = 1.0) was 5.4 m s−1 coming from the north. This led to a relatively high mean turbulent energy flux of 183 W m−2 (the turbulent flux at TMM with WWF = 1.0 was 182 W m−2). Figure 9 shows that areas of highest turbulent fluxes were on steep, windward (north-facing) slopes, while lower than average values can be detected on leeward slopes and in the valley bottoms. Overall, values range from 99 to 284 W m−2, with a standard deviation (S) of 8.8 W m−2, or 4.9% of the mean. Figure 10 shows the histogram of the distributed turbulent flux values.
The highest hourly values of turbulent fluxes were commonly observed during the afternoon hours as air temperatures and wind speeds typically reached a maximum. Studies have shown that wind speed over melting snowpacks tends to have a diurnal cycle, with maximum values in the early afternoon as the temperature inversion that often develops during the night weakens (Eaton and Wendler, 1982). Consequently, the highest absolute variabilities in turbulent fluxes were present in the basin during that time.
The hourly simulation results were summed to compute spatially variable daily values. Figure 11 illustrates the relative variability present in the daily values of turbulent fluxes by showing their standard deviation expressed as a percentage of the basinwide mean. The major factor determining daily relative variability over the study area was the amount of change in wind direction. Days with extremely variable wind directions (such as 22, 28, and 29 May 1999) produced low relative variabilities (Fig. 11). Higher spatial variability was found on days with constant wind directions. The highest relative variabilities were observed on days when sensible and latent heat fluxes changed directions during the course of the day. This, in conjunction with some change in wind directions, led to large differences in turbulent fluxes received in different parts of the basin (see Fig. 11; 14 May 1998 and 7 May 1999).
Daily values were summed to calculate accumulated spatially distributed turbulent fluxes for the entire modeling periods. Since northerly to northeasterly winds dominated the modeling period in 1999 (Fig. 8), the highest overall values can be seen on north- to northeast-facing slopes in the basin (Fig. 12). The lowest values are located on the steep (south-facing) leeward slopes, especially those leading toward the narrow river valleys and along the bottom of east-west–oriented valleys and gullies. Figure 13 shows the histogram of the computed values for the 105,000 grid cells of the DEM.
The standard deviation (S) for the accumulated turbulent fluxes in 1999 was found to be 2.3 MJ m−2 or about 3.4% of the area wide mean of 68 MJ m−2. To consider the importance of these small-scale variabilities for snowmelt, it is useful to convert the obtained turbulent fluxes at high and low wind speed locations into amounts of potential snowmelt. Assuming that all of the turbulent flux were used for snowmelt, locations at the opposite end of the 1S range would differ in their snowmelt amounts due to turbulent fluxes by 7 mm SWE. Areas within a range of 2S have as much as 9.2 MJ m−2 difference in their turbulent fluxes, leading to a difference in potential snowmelt of 28 mm SWE.
Model simulation for 1998 showed much the same results (not shown). Again, absolute variability was greatest during the afternoon hours, especially during the latter parts of the melt period, while relative daily variabilities again depended mainly on the amount of wind direction change. The distribution pattern of high and low amounts of accumulated turbulent fluxes was also very similar to 1999 (see Fig. 12) since both years had similar wind conditions. Overall, the standard deviation for the accumulated turbulent fluxes over the complete period was 6.0 MJ m−2, or 2.9% of the overall basinwide mean (205 MJ m−2). Expressed as potential snowmelt, this would translate to a difference of 18 mm SWE around the 1S range of values (72 mm for a range of 2S). Considering that the average end-of-winter SWE for the dominant open tundra areas of TVC typically fall between 50 and 120 mm, it becomes evident that the simulated differences in the turbulent fluxes and the resulting variable melt rates are a major contributor to the quick development of a patchy snow cover after the onset of melt.
Discussion
Spatial differences in turbulent fluxes have important implications for daily melt rates of snowpacks. Studies have shown that, after initiation of snowmelt at the top of the snowpack, the percolating meltwater forms two distinct wetting fronts with a faster moving finger front advancing ahead of a background front (Marsh and Woo, 1984). The input of spatially distributed snowmelt rates and snow cover depths is crucial for accurately simulating the advance of the two wetting fronts in different parts of the basin as was done by Marsh and Pomeroy (1996).
Turbulent fluxes also affect the timing of melt. Snow surface temperatures during the night, especially early in the melt season, are often lower than air temperatures. Consequently, turbulent fluxes, especially sensible heat fluxes, are directed toward the snow cover during those nights, counteracting the nighttime loss of radiative (long-wave) energy by the snow cover. This process was evident during both modeling periods in this study. As a result, snowpacks in areas of higher than average turbulent fluxes often experience a lower energy deficit during the night and warm up quicker during daylight hours, causing them to start to melt earlier in the day.
It is important for studies of snowmelt runoff from arctic catchments to correctly identify areas that have become snow-free, since these areas cease to contribute meltwater to runoff. Areas in the basin that have a low end-of-winter SWE in combination with increased amounts of melt energy will contribute earliest to snowmelt runoff and will subsequently become snow-free first. Variations in wind speed play a prominent role in both these processes. Studies in TVC and other arctic basins have shown that the snow cover is eroded over the winter, especially on slopes facing the predominant winter wind directions due to transport and sublimation losses during blowing snow events (Pomeroy and Li, 2000; Essery et al., 1999; Liston and Sturm, 1998; Pomeroy et al., 1997). The topography-influenced wind speed patterns lead to a positive feedback at high wind speed locations, combining eroded and therefore thinner end-of-winter snow covers with higher than average turbulent fluxes, provided the dominant wind direction does not change dramatically between the winter and spring periods. Such a pattern can be observed on north- and northwest-facing slopes in TVC.
The opposite scenario is simulated for locations of snow drifts, which may hold up to 33% of the entire TVC end-of-winter SWE (Marsh and Pomeroy, 1996). Above-average snow accumulations in TVC were observed for areas of low wind speeds and air flow separation such as the lee of hills, ridges, and uplands as well as for areas with higher shrub and forest vegetation and incised valleys (Pomeroy et al., 1997). These are also areas that receive below average contributions of turbulent fluxes throughout the melt period, leading to late-lying snow fields that persist well into the summer. Studies have illustrated the importance of such late-lying snow patches for the runoff regime of rivers in the Arctic (Marsh and Woo, 1981) and in more temperate environments (Luce et al., 1998). Late-lying snow drifts also have a major impact on runoff-producing mechanisms as infiltrating meltwater keeps the downslope water table high and close to the surface in the zone of high hydraulic conductivities (Quinton and Marsh, 1999).
As the snow cover becomes patchy, sensible energy is increasingly transferred laterally from the much warmer bare ground patches to the remaining snow patches (Liston, 1995; Marsh et al., 1997). The efficiency of this transfer is dependent on the size distribution of snow and bare ground patches but also on wind speed, increasing with higher wind speeds (Neumann and Marsh, 1998). This is another process where wind speed variations contribute to spatially variable melt rates.
The vast surface temperature differences between bare ground and snow covered surfaces, which can be as high as 42°C for arctic regions (Liston, 1995), cause the energy balances of snow and bare ground patches to be profoundly different. They should, therefore, be addressed separately in regional climate and numerical weather prediction models. These models often operate at scales that would treat the model area of this study (14 × 12 km) as one grid cell. Commonly, separate energy balances for snow covered and snow free areas are computed and grid cell averages from weighted means of the fractional areas of each land cover surface are subsequently calculated (Claussen, 1991). However, this “tile model” approach does not consider the small-scale turbulent flux differences over snow covered areas presented in this study. Furthermore, the approach is not well suited for calculating average turbulent fluxes over a larger area with pronounced surface inhomogeneities (i.e., snow/no snow) since it does not consider any interaction between the tiles. Therefore, it cannot account for local advection which plays an important role in the overall energy flux of the composite landscape (Marsh, 1999; Lhomme et al., 1994).
Conclusions
This study examined the topographic influence on surface windflow and the resulting impact on the small-scale variability of sensible and latent heat fluxes during snowmelt of an arctic catchment. A detailed field experiment showed that there are pronounced differences in wind speeds introduced by the topography, even in the relatively low relief study area. Wind speed differences of up to 39% were found between high wind speed locations (the top and the upper windward slopes of a hill) and low wind speed areas on lower leeward hillslopes. Spatially variable surface wind fields could be modeled fairly accurately with a relatively simple windflow model and a high-resolution DEM. The model was validated with data from a field experiment and four meteorological stations within the basin.
Distributed hourly sensible and latent heat fluxes for 33 days during the springs of 1998 and 1999 were computed using a bulk aerodynamic approach and modeled distributed wind fields. The results show that the small-scale turbulent flux variability contributed greatly to spatially variable melt rates. Furthermore, variable wind speeds throughout winter and spring are a major factor in the development of a patchy snow cover as shallow end-of-winter snowpacks are combined with above-average turbulent energy fluxes in areas of increased wind speeds. Snow drifts develop in low wind speed areas over the winter and persist well into the summer partially due to decreased amounts of turbulent fluxes received during melt.
The study showed that turbulent fluxes play an important role in determining the location of snow covered areas and the timing and volume of meltwater contributions of these areas to runoff. The hydrologic and energetic important steeper slopes of the basin differed in their turbulent fluxes by more than 9 MJ m−2, or about 28 mm of potential snowmelt. Hydrological, snowmelt-runoff models should therefore include spatially distributed turbulent fluxes, especially if they operate on the concept of varying runoff source areas. Wind directions and locations of above- and below-average turbulent fluxes did not differ much between the two years. The relative variability of the studied fluxes was also fairly constant. This indicates that it might be possible to determine overall variabilities of turbulent fluxes in open, relatively flat landscapes from terrain data such as frequency distribution of slopes and aspects, which can be easily obtained from DEMs. Mean turbulent fluxes and their variability within an area could then be computed from point measurements of air and surface temperature, humidity, wind speed, and wind direction without having to run a detailed small-scale model. However, more work is needed in different landscapes and different terrain types to test this concept and to develop such statistical relationships between terrain properties and spatial variability of turbulent fluxes.
Acknowledgments
The authors would like to thank Cuyler Onclin, Natasha Neumann, Mark Russell, Ross MacKay, Bruce Davison, and Steve McCartney for their essential help in collecting and interpreting field data. Field support was provided by the Aurora Research Institute. This study received funding and logistical support from the National Water Research Institute, the Canadian GEWEX/MAGS program, and the Polar Continental Shelf Project. The authors are grateful for the comments of two anonymous reviewers that greatly improved this manuscript.
References Cited
FIGURE 1.
Digital elevation model of Trail Valley Creek and surrounding areas. Also shown are the locations of the four meteorological stations operating in the basin during the spring of 1999. Trail Valley Main Met (TMM) and Trail Valley Upper Plateau (TUP) are long-term stations, while Trail Valley Deep Snow (TDS) and Trail Valley Bare Ground (TBG) are temporary stations located over deep and shallow snow sites, respectively. Box A denotes the field experiment shown in Figure 4
FIGURE 2.
Wind weighting factors for the digital elevation model for a wind direction of 45°. Box A indicates location of the anemometer experiment in 2002 (see Fig. 4)
FIGURE 3.
Calculated vs. observed snow surface temperatures at TDS in Trail Valley Creek during spring of 1996
FIGURE 4.
Instrumented hill with anemometer locations. Also shown are modeled wind weighting factors for a wind direction of 45° (Topographic contour interval = 10 m)
FIGURE 5.
Modeled vs. observed wind speeds at (a) Upper Plateau met station (TUP) for spring of 1999, and (b) station 5 located on south-facing hillslope for spring period of 2002
FIGURE 6.
Average daily meteorological conditions in Trail Valley Creek for 5 May to 6 June, 1998 and 1999: (a) daily temperature (Ta), (b) relative humidity (RH), and (c) windspeed
FIGURE 7.
Average daily accumulated sensible and latent heat fluxes in MJ m−2 for (a) 5 May to 6 June 1998 and (b) 5 May to 6 June 1999
FIGURE 11.
Standard deviations of accumulated daily turbulent fluxes throughout modeling area expressed as percentage of area-wide mean
